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Abstract. We have developed an analytic model to describe coupling of plasma 
and neutral fluids in the partially ionized heliosphere plasma medium. The sources 
employed in our analytic model are based on a K-distribution as opposed to the 
Maxwellian distribution function. Our model uses the K-distribution to analytically 
model the energetic neutral atoms that result in the heliosphere partially ionized 
plasma from charge exchange with the protons and subsequently produce a long tail 
which is otherwise not describable by the Maxwellian distribution. We present our 
analytic formulation and describe major differences in the sources emerging from 
these two distinct distributions. 



1. Introduction 

With the Voyager spacecraft now in the heliosheath (see Fig 1), the in situ character 
of the solar wind plasma can be explored. Surprisingly, the supersonic solar wind 
plasma, probed by the ACE/WIND/Cluster spacecrafts near 1 AU (Astronomical 
Unit), depicts an entirely different character when contrasted with the Voyager I 
and 2 observations in the heliosheath region (typically beyond 84 AU) (Goldstein 
et al 1995, Burlaga et al 2005, 2006, 2008, 2009; Stone et al 2005; Decker et al 2005; 
Richardson et al 2008, Zank 1999). Figure 1 shows an idealized cartoon reflecting 
our current understanding based on theory, simulations and modeling together with 
observations. Little is known about the physical processes that govern the intricate 
multiscale (associated with waves, structures, turbulence) interactions outlined in 
Fig 1. The supersonic solar wind (SW) plasma interacts with local interstellar 
medium (LISM) neutral hydrogen (H) gas through charge exchange leading to the 
creation of energetic pick up ions (PUI). The SW is decelerated, compressed and 
heated at a shock, the termination shock (TS), across which it develops small scale 
turbulence (Shukla 1978, Shaikh & Zank 2008, 2010, Shaikh 2010, Shaikh et al 2006, 
Mendonca & Shukla 2007). In the heliosheath region, the nonlinear structures, such 
as magnetic hole and humps are found (Burlaga et al 2008, 2009). The SW pro- 
tons continue to interact with neutrals via charge exchange to produce significant 
number of pick up ions. Both in the supersonic and subsonic SW (at least in the 
outer heliosphere) the pressure associated with the PUIs exceeds that of the solar 
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wind protons (Burlaga et al 2009). Both Voyager 1 and 2 are reporting a number of 
puzzling observations that were not anticipated by existing analytic or simulation 
models. An intriguing example is that of magnetic field distribution. The latter 
is lognormal in supersonic solar wind, whereas it exhibits a Gaussian distribution 
in subsonic heliosheath. Surprisingly, Voyager 2 indicates that the magnetic field 
distribution is lognormal in the subsonic heliosheath plasma. The source of this 
apparent discrepancy in the magnetic field distributions reported by Voyagers 1 & 
2 in the heliosheath is not known. Another example is that of plasma in heliosheath 
which is comprcisscd, turbulent and it is an admixture of waves, fluctuations and 
magnetic structures (magnetic hole/hump, see section 2 for details) (Burlaga 2006, 
2009). The effect of PUIs on the formation and evolution of nonlinear magnetic 
structures, waves and fluctuations in outer heliosphere and the heliosheath plasma 
are an open question. These issues continue to pose severe challenges to our under- 
standing of the heliosheath plasma. 

Although there exists wealth of in situ measurements by the Voyager space- 
crafts, they do not provide much information about the global structure of the 
heliosphere interactions. For instance, the coupling of plasma protons with the in- 
terstellar neutral atoms has traditionally been done through MaxwcUian sources. 
However, careful studies have revealed that the distribution of hydrogen neutral 
(after charge exchanging, they turn into energetic neutral atoms, ENA) does not 
exactly follow a Maxwellian functions. Recently, Prested et al. (2008) used a k- 
distribution for the ENA parent population to obtain ENA maps. The advantage 
of using this distribution, as opposed to a Maxwellian, is that it has a power-law 
tail, and is therefore capable of producing ENAs at suprathermal energies. Now 
there has been an increasing consensus that the plasma and neutral fluids follow 
nearly kappa distribution (Heerikhuisen et al 2008). 

A realistic modeling of the heliosheath plasma, one that includes a self-consistent 
treatment of the PUIs, is therefore critically important and essential to our un- 
derstanding of the highly variable heliosphere plasma. The central them of this 
paper is therefore to model complex coupling between plasma and neutral fluids 
via K-distribution as oppose to the Maxwellian distribution. Note here that the 
K distribution modifies the charge exchange interactions in fluid equations. The 
kappa-distribution emphasizes charge exchange by high temperature protons. We 
will investigate the efl^ects of the K-distribution in heliosphere plasma turbulence 
for single fluid plasma-neutral coupled turbulence models. 

In section 2, we describe K-distribution for neutral and plasma distribution and 
derive sources for the complex coupling interactions between the two distinct fluids. 
Section 3 describes complete source terms for the coupling interactions. Finally, a 
summary is presented in section 5. 

2. Plasma neutral coupling via ^-distribution source 

The charge exchange terms can be obtained from the Boltzmann transport equa- 
tion that describes the evolution of a neutral distribution function /(r,v,<) in a 
six-dimensional phase space defined respectively by position and velocity vectors 
{x,Vx-iVy,Vz) at each time t. Here we follow Pauls et al. (1995) in computing the 
charge exchange terms, based on K-distribution functions, from various moments 
of the Boltzmann equation. The Boltzmann eqiiation for the neutral distribution 
contains a source term proportional to the proton distribution function fp and a 
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Heliosphere: region filled with solar wind plasma and bounded by HP. 



HT 




HT: Heliotail 
HP: Heliopause 
IHS: Inner Heliosheath 
OHS: Outer Heliosheath 
BS: Bow Shock 
TS: Termination Shock 
HO: Magnetic hole 
HU: IMagnetic hump 



Structures 



Waves 



Turbulence 

Figure 1. Schematic overview of different regions in the global heliosphere. The solar 
wind emanating from the Sun propagates outward and interacts with partially ionized 
interstellar gas predominantly via charge exchange, and creates pick up ions (PUIs). At 
the termination shock (TS), the supersonic SW decelerated, heated, compressed becoming 
subsonic, in the heliosheath, and again interacts with interstellar neutrals via charge ex- 
change before it reaches heliopause (HP). The subsonic SW flows down into the heliotail. 
Magnetic structures such as magnetic holes/humps are observed in heliosheath plasma. 
During its journey from the Sun to the IfP, the solar wind plasma develops multitude of 
length and time scales that interact with the partially ionized interstellar gas, TS, and 
nonlinear structures develop in a complex manner. 



loss term proportional to the neutral distribution function /„ . 

4/p(r, V, t) + V V/p(r, V, i) + - • V,/p(r, v, t) 
ot m 

/„(r,V,i) j fpir,Wp,t)\-Vp-w\(Jex{Vrel)d?^p- 

fp{r,v,t) / /„(r,v„,t) I {Vrel)d 



(2.1) 

where, Uex{vrei) is the charge exchange cross section. The charge exchange param- 
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of the neutrals and the protons through aex = [(2.1 — 0.092 ln(urei))10~'' cm]^ (Fite 
et al 1962). This cross-section is valid as long as energy does not exceed leV, which 
usually is the case in the inner/outer heliosphere. Beyond leV energy, this cross- 
section yields a higher neutral density. This issue is not applicable to our model 
and hence we will not consider it here. The density, momentum, and energy of the 
thermally equilibrated MaxwcUian proton and neutral fluids can be computed from 
Eq. (2.1) by using the zeroth, first and second moments / f^d^^, J m^f^d'^S, and 
/ m^^/2/j(i'^^ respectively, where ^ = Up or w„. Since charge exchange conserves 
the density of the proton and neutral fluids, there are no sources in the correspond- 
ing continuity equations. We, therefore, need not compute the zeroth moment of 
the distribution function. Computing directly the first moment from Eq. (2.1), we 
obtain the neutral fluid momentum equation. 

A similar evolution equation can be written for the neutral distribution function 
/n(i")Vn,t). We consider the case where both fp{r,Vp,t) and /„(r, Vn,t) are given 
by a K distribution of the following type: 



/p(r,Vp) 



Up T{k+1) 



1 + 



(Vp 



Up)^ 



1 



KV^ 



(2.2) 



/n(r, v„) 



n„ r(K+l) 



T^^Vj,^ KsF (k — I) 

First we evaluate the following integral: 



(v„ - U„)2 



(2.3) 



Pp{r,V,t) = (Jex{Vrel) j fp{r,Vp,t) |Vp - v| rf^Vpi (2.4) 

where (7ex{vrei) is taken out of the integral, as it varies slowly with respect to {vrei)- 
The integral (2.4) is fully written as 



Tl 



1 



with 



We write, 

Vp - Up = (vp - v) - (Up - v) 

and define new variables as 



|Vp-v|d3vp; (2.5) 



V = (vp- v)/y'Ki;2^; x = (Up - v)/.^^^^; 
with the new variables, the integral in Eq (2.5) becomes. 



Pp{r,V,t) = aex{Vrel)- 



/oo 
[l + (V-x)2] 
-cx> 



Vd-'V: 



(2.6) 
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where we have used, 

|vp - v| = yfi^V; d^Vp = {nv'^^f'^d^Y 
and the constant before the integral in Eq (2.6) is, 



T{n + 1) HpVT^ T{n + 1) 



We now proceed to evaluate the integral Eq (2.6). In spherical coordinate, 

d^V = V^dvsme de d(j), 

where 6 is the angle between V and x, after performing the ^ integration, with 
fj. = cos 6 

/ = 27rJ^ V^dV J dn{l + V^ -2Vxn + x^)-^''+'^^ (2.7) 
This integration becomes, 

" pX pX pOO 

1 = 2 / z'^il + z'^ydz + x'^ {l + z'^)-''dz + 2x z{l + z'^ydz 

Jo Jo Jx 

with X = (Up — v)/ ^Kv^ . We now proceed to determine the explicit values of the 

above definite integrals. The first two integrals are given in terms of the hypergeo- 
metric functions, 2F1, which are 

£ z^(l + zY'^dz=^ ,F,(^l,K;l;-x'y 



x^ {l + z^)-''dz = x 2Fl(^^,K■,^■,-x' 



where the Hypcrgcomctric function 2F1 (a, 6; c; z) (with a, b, c are constant numbers 
and z is the variable) is expressed as a power series in z: 

, , ab z , a{a+l)b{b+l) z^ 

2Fi{a,b; c; z) = 1 + — - + — — — 

c 1! c(c + 1) 2! 

, a(a+l)(a + 2)6(6+l)(6 + 2) z3 

+ c(c+l)(c+2) 3!+- ^'-^^ 

Using Kummer identity for hypergeometric functions, 

2F1 (a, b; c; z) = 2F1 (6, a; c; -j) = (1 - z)-'>2Fi [b, c-a;c; z/{z - 1)] (2.8a) 



2F, [l^; l,-xA = {iWrSF^ ) = (l+x^)-«2Fi I 1, ^ 



(2.86) 

Similarly, 

2F1 I -xA = (1 + x2)-«2Fi (1, «; ^; (2.8c) 
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The last integral can be evaluated easily 



f 

J X 



z{i + z'^y^dz 



(1 + 



2(-K+l) 



-K+l 



2(fi;- 1) 



(k must be > 1) 



Collecting the above terms, the integral /3p(r, v) is, 



TTK X 

,2 



r(«-i) 



+ — 2F1 



2' l + a;2 



+ 



1+^2 



K- 1 



; a;=|Up-v|/y^ 



(2.9) 

An approximate value of the above expression in the two limits ^/k a; <C 1 and 
a; 1 can be obtained as follows: 



a; «; 1 : /? = 



2npaexVTp T{k + 1) 



TTK r(«-i) 



1 ^ (Up - v) 



K-1 



^1 



.X > 1 : 



I ^VT.T^jn + 1) 

7rK(K- l)r2 (k- £ 



(2.10) 
(Up-v)2 (2.11) 



Note that in an asymptotic limit, the Gamma functions for large argument is 

lim r{K + a) ~ k'T(k) 



lim k" — — rr 1, 



K^-oo T{k + b) 



3. Complete expressions for the source terms 

To find the source terms, wc take the moments of Eq (2.1) by multipling both sides 
with various powers of the velocity v. The zeroth order moment would contribute 
to the source term of the mass continuity equation, the first order moment would 
contribute to the source term of the momentum equation, the second order moment 
would contribute to the source term of the energy equation. The moments for the 
left hand terms of Eq (2.1) are well known, so we shall show the moments of the 
right hand terms, using kappa distribution for fp and /„. We derive full expression 
for the integrals on the r.h.s of Eq (2.1) by using the complete expression for the 
integral for /3p(r, v, t) given by Eq (2.9) without any approximation. 

Since the charge exchange process conserves the proton and neutral density, there 
will be no source term for the mass continuity term, so we need not calculate the 
zeroth moment. Hence we start with the first moment of the right hand side of Eq 
(2.1). With 

= /n(r,v-U„), 

and /3p(r, v, t) = f3p{r, x, t) = f3p{T, Up - v) 

The production term for Momentum transport, is, 



/n(r,v) = 



r(K + 1) 



(v - U„)^ 
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/■°° ^ 

Qmp = d w/„(r,v - U„)^p(r,Up - v) 

Jo 

1*00 

= / d3fU„/„(r,v-U„)/3p[r,(Up-U„)-(v-U„)] 
Jo 

+ / d?v{w - U„)/„(r, V - U„)^p[r, (Up - U„) - (v - U„)], 

Jo 

foo 1*00 

= U„/ d3zi/„(r,u)/3p[r,(AU-u)]+ / d3Mu/n(r, u)/3p[r, (AU - u)]; 
Jo Jo 



(3.1) 

where we used, u = v — Un; AU = Up — Un- 
The production term for Energy transport, is, 

/>oo 

Qep= / d3i;|v|2/n(r,v-U„)^p(r,Up-v) 
Jo 

/•OO 

= / rf=^^;|v-U„|Vn(r,v-U„)/3p(r,Up-v) 

/>oo 

+2U„- / dWn(r,v-U„)/3p(r,Up-v) 

/■OO 

-C/^ / d^t;/„(r, V - U„)/3p(r, Up - v) 
Jo 

= / d3^|v-U„|2/„(r,v-U„)/3p(r,Up-v) 
Jo 

POO 

+2U„ • / - U„)/„(r, V - U„)/3p(r, Up - v) 

/>oo 

+2C/2 / (i3i;/„(r, V - U„)^p(r, Up - v) 

Jo 

POO 

-Ul / d^vUv, V - U„)/3p(r, Up - v) (3.2) 
Jo 

As before, we introduce the variables u = v — U„; AU = Up — Un and write 
/3p(r, Up - v) = /3p(r, (Up - Un) - (v - Un)) = /3p(r, AU - u). Expression (3.2) 
can be written as 

POO poo 

Qep= / dW/„(r,u)/3p(r,AU-u)+2Un- / d\ufn{v,u)l3p{r, AV - u) 
Jo Jo 

POO 

+Ul / d^'ufniv, u)/3p(r, AU - u) (3.3) 
Jo 



4. Summary and conclusion 

In summary, our major results are Eq (3.1) & Eq (3.2). A tentative comparison of 
the sources based on Maxwellian and kappa distribution functions is shown in Fig 

(2). It is evident from this figure that the sources based on a Maxwellian distribution 
function falls off sharply and without any tail region. This therefore excludes the 
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Figure 2. Comparision between the kappa and Maxwellian distribution functions for the 
sources. Clearly, the tail of the distribution for k function is long and wide as opposed to 
the Maxwellian distribution. It is because of this feature, we have computed sources based 
on a K-distribution function. 

energertic component of the neutral atoms and hence is inappropriate for typical 
ENAs. By contrast, the sources based on a kappa distribution function depicts a 
well-behaved tail distribution that represents ENA distribution. 

Our previous work in Shaikh & Zank (2008) has shown that charge exchange 
modes modify the helioshperic turbulence cascades dramatically by enhancing non- 
linear interaction time-scales on large scales. Thus the coupled plasma system 
evolves differently than the uncoupled system where large-scale turbulent fluctua- 
tions are strongly correlated with charge-exchange modes and they efficiently be- 
have as driven (by charge exchange) energy containing modes of helioshperic tur- 
bulence. By contrast, small scale turbulent fluctuations are unaffected by charge 
exchange modes which evolve like the uncoupled system as the latter becomes less 
important near the larger k part of the helioshperic turbulent spectrum. The neutral 
fluid, under the action of charge exchange, tends to enhance the cascade rates by 
isotropizing the helioshperic plasma turbulence on a relatively long time scale. This 
tends to modify the characteristics of helioshperic plasma turbulence which can 
be significantly different from the Kolmogorov phenomenology of fully developed 
turbulence. It remains to be seen how these modified sources influence nonlinear 
turbulent properties of the small scale helioshperic plasma fluctuations. 
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